home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / obsolete / goodfit.pro < prev    next >
Text File  |  1997-07-08  |  4KB  |  174 lines
  1. ; $Id: goodfit.pro,v 1.2 1997/01/15 04:02:19 ali Exp $
  2. ;
  3. ;  Copyright (c) 1991-1997, Research Systems Inc.  All rights
  4. ;  reserved. Unauthorized reproduction prohibited.
  5.  
  6.  
  7.  
  8. Function Stdev1,D,Freq
  9. ; Stdev returns the sample standard deviation, where
  10. ; D contains the sampled values and Freq their frequency.
  11. FT=Total(Freq)
  12. D2=Total(D*Freq) / FT
  13. Return, SQRT(Total(freq*(D - D2)^2)/(FT-1))
  14. END    ;Function Stdev
  15.  
  16.  
  17.  
  18. Pro goodfit, Freq,A, B,ChiSqr, Prob,DF,Distr=D
  19. ;+
  20. ; NAME:
  21. ;       GOODFIT
  22. ;
  23. ; PURPOSE:
  24. ;    Test that a set of data has a given distibution.  User can select
  25. ;    built-in distribution through the keyword DISTR or supply own 
  26. ;    expected values.
  27. ;
  28. ; CATEGORY:
  29. ;    Statistics.
  30. ;
  31. ; CALLING SEQUENCE:
  32. ;    GOODFIT, Freq, [A, B, ChiSqr, Prob, DF, DISTR=D]
  33. ;
  34. ; INPUTS: 
  35. ;     Freq:    Vector of value or range frequencies.
  36. ;          
  37. ; OPTIONAL INPUTS:
  38. ;     A:    Vector of observed values or left endpoints for interval data.
  39. ;        Should have same length as F.  A must be present if user
  40. ;        selects built-in distributions.                
  41. ;
  42. ;     B:    Vector of right hand endpoints.
  43. ;                            
  44. ;
  45. ; OPTIONAL OUTPUT PARAMETERS:
  46. ;    ChiSqr:    The chi square statistic to test for the goodness of fit of
  47. ;        the data to the distribution.
  48. ;
  49. ;    Prob:    The probability of ChiSqr or somethinG larger from a chi 
  50. ;        square distribution.
  51. ;  
  52. ;    DF:    Degrees of freedom. User must specify chi
  53. ;        square degrees of freedom for user supplied
  54. ;        distribution. Compute DF to be s-1-t,
  55. ;        where s = N_Elements(A) and t is the number
  56. ;        of parameters that must be estimated to
  57. ;        derive expected frequencies.
  58. ;
  59. ; KEYWORDS: 
  60. ;    Distr:    the type of distribution for which the data is to be tested.
  61. ;        If Distr = "G", then GOODFIT tests for a Gaussian distribution.
  62. ;        To supply distribution, set Distr to the vector of expected
  63. ;        frequencies, say, EF.  EF(i) = the frequency expected for the
  64. ;        observation with observed frequency Freq(i).
  65. ;
  66. ; RESTRICTIONS:
  67. ;    None.
  68. ;
  69. ; COMMON BLOCKS:
  70. ;    None.
  71. ;
  72. ; SIDE EFFECT:
  73. ;    None.
  74. ;
  75. ; PROCEDURE:
  76. ;    Compute expected frequencies EXP(i).
  77. ;    Chi_Sqr = SUM((Freq(i)-ExpV(i))^2/ExpV(i)^2)
  78. ;    has the chi sqr distribution.
  79. ;-
  80. On_Error,2
  81. SF= SIZE(Freq)
  82.  
  83.  
  84. if(KEYWORD_SET(D) eq 0) THEN BEGIN
  85.    print,"Goodfit- key word Distr must be set"
  86.    Return
  87.  ENDIF ELSE if (N_ELEMENTS(D) gt 1) THEN BEGIN
  88.               Expv = D
  89.               Distr = "U"
  90.              ENDIF $
  91.            ELSE Distr = D
  92.  
  93.  if Distr NE "U" THEN  BEGIN
  94.    SA= Size(A)
  95.    if((SA(0) NE 1) OR (SF(0) NE 1) OR (SF(1) NE SA(1))) THEN BEGIN
  96.      print, "GoodFit- incompatible vectors"
  97.     Return
  98.    ENDIF
  99.  ENDIF
  100.  
  101.  
  102.  
  103. Case Distr of
  104.  
  105.  
  106. "G": BEGIN
  107.  
  108.      if(N_Elements(B) EQ 0) THEN BEGIN
  109.      print,"GoodFit- interval right endpoints undefined"
  110.      RETURN
  111.      ENDIF
  112.  
  113.  
  114.      D=((A+B)/2.0)
  115.      FT=Total(Freq)
  116.  
  117.      if(FT LE 1) THEN BEGIN
  118.         print,"GoodFit-must have more than 1 observed value" 
  119.         return
  120.      END
  121.    
  122.      mean=Total(D*Freq)/FT
  123.      sigma=StDev1(D,Freq)
  124.      A1=(A-mean)/sigma
  125.      B1 = (B-mean)/sigma
  126.      ExpV=FltArr(SA(1))
  127.      for i=0,SA(1)-1 DO                       $
  128.        Expv(i) = (Gaussint(B1(i))-Gaussint(A1(i)))*FT
  129.    
  130.      DF=SA(1)-2
  131.      END
  132.  
  133.  "U": BEGIN
  134.  
  135.     if(N_Elements(Distr) EQ 0) THEN BEGIN
  136.         print,"GoodFit- expected values undefined"
  137.         RETURN
  138.      ENDIF
  139.  
  140.      if (N_Elements(DF) EQ 0 ) THEN BEGIN
  141.       print, "goodfit- must specify degrees of freedom"
  142.       Return
  143.      ENDIF
  144.  
  145.      SE=Size(EXPV)
  146.      if(SE(0) NE 1) OR (SE(1) NE SF(1)) THEN BEGIN
  147.        print, "GoodFit- incompatible arrays"
  148.        Return
  149.      ENDIF     
  150.      END
  151.  
  152.   
  153.  ELSE: BEGIN
  154.        print,"Unknown value for Key word Distr"
  155.        Return
  156.        END
  157.  
  158.  ENDCASE
  159.  
  160.  SS=float(Freq-Expv)
  161.  ChiSqr = Total(SS*SS/Expv)
  162.  Prob = 1 - chi_sqr1(ChiSqr,DF)
  163.  print, "ChiSqr = ",ChiSqr, " with  ",DF,             $
  164.    " degrees of freedom"
  165.  print, "Probability of ChiSqr=",Prob
  166.  
  167.  Return
  168.  END
  169.  
  170.    
  171.  
  172.  
  173.  
  174.